cbPalette <- c("#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7")

plot.performance <- function(alloc, pred, se, test, model.name) {
	# Plots the performance of each method.
	# Args:
	#   alloc: the data frame of portfolio weight allocation
	#   pred:  the data frame of predicted returns of portfolios during the test period
	#   se:    the data frame of the predicted se during the test period
	#   test:  the test data
	png(paste("../results/", model.name, "_predict.png", sep=""))
	par(mfrow=c(3,2))
	
	# Plot the prediction region.
	#smlo
	smlo.pred <- ts(pred$smlo, freq=12, start=c(1997,1))
	smlo.se <- ts(se$smlo, freq=12, start=c(1997,1))
	smlo.truth <- ts(test$smlo_vwret,  freq=12, start=c(1997,1))
	ts.plot(smlo.pred, smlo.pred+smlo.se, smlo.pred-smlo.se, gpars=list(lty=c(1,2,2), ylab="smlo"))
	points(smlo.truth, col="red")
	correct = sum((smlo.truth <= smlo.pred + smlo.se) & (smlo.truth >= smlo.pred - smlo.se))
	wrong = length(smlo.truth) - correct
	smlo.prec = correct / length(smlo.truth)
	
	#smme
	smme.pred <- ts(pred$smme, freq=12, start=c(1997,1))
	smme.se <- ts(se$smme, freq=12, start=c(1997,1))
	smme.truth <- ts(test$smme_vwret,  freq=12, start=c(1997,1))
	ts.plot(smme.pred, smme.pred+smme.se, smme.pred-smme.se, gpars=list(lty=c(1,2,2), ylab="smme"))
	points(smme.truth, col="red")
	correct = sum((smme.truth <= smme.pred + smme.se) & (smme.truth >= smme.pred - smme.se))
	wrong = length(smme.truth) - correct
	smme.prec = correct / length(smme.truth)
	
	#smhi
	smhi.pred <- ts(pred$smhi, freq=12, start=c(1997,1))
	smhi.se <- ts(se$smhi, freq=12, start=c(1997,1))
	smhi.truth <- ts(test$smhi_vwret,  freq=12, start=c(1997,1))
	ts.plot(smhi.pred, smhi.pred+smhi.se, smhi.pred-smhi.se, gpars=list(lty=c(1,2,2), ylab="smhi"))
	points(smhi.truth, col="red")
	correct = sum((smhi.truth <= smhi.pred + smhi.se) & (smhi.truth >= smhi.pred - smhi.se))
	wrong = length(smhi.truth) - correct
	smhi.prec = correct / length(smhi.truth)
	
	#bilo
	bilo.pred <- ts(pred$bilo, freq=12, start=c(1997,1))
	bilo.se <- ts(se$bilo, freq=12, start=c(1997,1))
	bilo.truth <- ts(test$bilo_vwret,  freq=12, start=c(1997,1))
	ts.plot(bilo.pred, bilo.pred+bilo.se, bilo.pred-bilo.se, gpars=list(lty=c(1,2,2), ylab="bilo"))
	points(bilo.truth, col="red")
	correct = sum((bilo.truth <= bilo.pred + bilo.se) & (bilo.truth >= bilo.pred - bilo.se))
	wrong = length(bilo.truth) - correct
	bilo.prec = correct / length(bilo.truth)
	
	#bime
	bime.pred <- ts(pred$bime, freq=12, start=c(1997,1))
	bime.se <- ts(se$bime, freq=12, start=c(1997,1))
	bime.truth <- ts(test$bime_vwret, freq=12, start=c(1997,1))
	ts.plot(bime.pred, bime.pred+bime.se, bime.pred-bime.se, gpars=list(lty=c(1,2,2), ylab="bime"))
	points(bime.truth, col="red")
	correct = sum((bime.truth <= bime.pred + bime.se) & (bime.truth >= bime.pred - bime.se))
	wrong = length(bime.truth) - correct
	bime.prec = correct / length(bime.truth)
	
	#bihi
	bihi.pred <- ts(pred$bihi, freq=12, start=c(1997,1))
	bihi.se <- ts(se$bihi, freq=12, start=c(1997,1))
	bihi.truth <- ts(test$bihi_vwret,  freq=12, start=c(1997,1))
	ts.plot(bihi.pred, bihi.pred+bihi.se, bihi.pred-bihi.se, gpars=list(lty=c(1,2,2), ylab="bihi"))
	points(bihi.truth, col="red")
	correct = sum((bihi.truth <= bihi.pred + bihi.se) & (bihi.truth >= bihi.pred - bihi.se))
	wrong = length(bihi.truth) - correct
	bihi.prec = correct / length(bihi.truth)	
	dev.off()
	
	# Plot the prediction precision
	png(paste("../results/", model.name, "_precision.png", sep="")) 
	barplot(c(smlo.prec, smme.prec, smhi.prec, bilo.prec, bime.prec, bihi.prec), ylab="Precision of prediction",
		names.arg=c("smlo", "smme", "smhi", "bilo", "bime", "bihi"), ylim=c(0, 1)) 
	dev.off()
	
	# Plot the transition map of the weight allocation
	plot.abs.transition.map(alloc[,2:7], paste("../results/", model.name, "_transition.png", sep=""))
}

plot.acf <- function() {
	png("../results/vwret_acf.png")
	par(mfrow=c(3,2))
	acf(FF6$smlo_vwret)
	acf(FF6$smme_vwret)
	acf(FF6$smhi_vwret)
	acf(FF6$bilo_vwret)
	acf(FF6$bime_vwret)
	acf(FF6$bihi_vwret)
	dev.off()
}

plot.returns <- function(FF6, filename) {
	# Plots the return series in the data in different charts.
	# Args:
	#   data: the data frame containing all the returns for smlo, smme, smhi, bilo, bime, bihi
	#   filename: the png filename to save the charts to.
	require(ggplot2)
	png(filename)
	q1 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="red", aes(y=smlo_vwret)) +
		xlab("date") + ylab("portfolio returns for smlo") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	q2 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="orange", aes(y=smme_vwret)) +
		xlab("date") + ylab("portfolio returns for smme") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	q3 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="yellow", aes(y=smhi_vwret)) +
		xlab("date") + ylab("portfolio returns for smhi") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	q4 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="green", aes(y=bilo_vwret)) +
		xlab("date") + ylab("portfolio returns for bilo") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	q5 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="blue", aes(y=bime_vwret)) +
		xlab("date") + ylab("portfolio returns for bime") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	q6 <- ggplot(FF6, aes(x=as.Date(as.character(date), format="%Y%m%d"))) +
		geom_line(colour="purple", aes(y=smlo_vwret)) +
		xlab("date") + ylab("portfolio returns for bihi") +
		scale_x_date() + scale_y_continuous(limits=c(-0.3,0.3))
	multiplot(q1, q2, q3, q4, q5, q6, cols=2)
	dev.off()
}

plot.transition.map <- function(data, filename) {
	# Plot data in a transition map, where x is the time, y is stacked weights
	# Args:
	#   data: the data frame containing all the weights for smlo, smme, smhi, bilo, bime, bihi
	#   filename: the png filename to save the charts into.
	png(filename)
	barplot(t(as.matrix(data)), col=c("#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00"), legend = c("smlo", "smme", "smhi", "bilo", "bime", "bihi"), ylab="Portfolio Weights", xlab="Time")
	dev.off()	
}

plot.abs.transition.map <- function(data, filename) {
  # Plot data in a transition map, where x is the time, y is stacked weights
  # Args:
  #   data: the data frame containing all the weights for smlo, smme, smhi, bilo, bime, bihi
  #   filename: the png filename to save the charts into.
  png(filename)
  barplot(t(as.matrix(abs(data))), col=c("#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00"), legend = c("smlo", "smme", "smhi", "bilo", "bime", "bihi"), ylab="Portfolio Weights", xlab="Time")
  dev.off()	
}

compute.returns <- function(alloc, test, are.log) {
	# Compute the return for each time period based on the allocation weight and the real returns:
	# Args:
	#   alloc: the data frame of portfolio weight allocation
	#   test: the test data set with real return data
  #   are.log: logical parameter whether returns are log trasnformed
	# Returns:
	#   a time series with cumulative returns
	wt <- alloc[,2:7]
	vwret <- data.frame(smlo=test$smlo_vwret, smme=test$smme_vwret, smhi=test$smhi_vwret, bilo=test$bilo_vwret, bime=test$bime_vwret, bihi=test$bihi_vwret)
	r <- rowSums(wt * vwret) 
  	# @conxing. I dont think you can sum returns unless they are log. you want cumprod.
  	# see my get.cumulative.returns() below as an example.
  	if(isTRUE(are.log)){
  		r.cum <- cumsum(r)
  	} else{
    	r.cum <- cumprod(r +1) -1
  	}
	return(ts(r.cum, freq=12, start=c(1997,1)))
}

get.cumulative.returns <- function(returns.data, are.log)
{
  # are.log is a logical expression. returns are in same format as supplied
  if(isTRUE(are.log)){
    cumulative.returns <- (cumsum(returns.data))
  } else{
    cumulative.returns <- (cumprod(returns.data + 1) -1) # this is for non-log returns
  }
  return(cumulative.returns)
}

# please review this since the results are not quite matching my expectations for calculations by hand
library("tseries")
get.sharpe.ratio <- function(cumulative.returns, rf.rate=0) {
	sharpe <- rep(NA, length(cumulative.returns))
	for (i in 1:length(sharpe)) {
		sharpe[i] <- monthly.sharpe.ratio(cumulative.returns[1:i], rf.rate) 
	}
	return(sharpe)
}

monthly.sharpe.ratio <- function(cumulative.returns, rf.rate=0){
  # returns must be converted to cumulative
  my.sharpe.ratio <- sharpe(ts(cumulative.returns, frequency=12), r=rf.rate, scale = sqrt(12))
  return(my.sharpe.ratio)
}

